Overview of analysis

This analysis examines if overall community composition differs based on self-reported ethnicity and what specific taxa are differentially abundant. MetaPhlAn (v3.0.7) was used to profile the taxonomic composition of metagenomes and calculate unweighted and weighted UniFrac distances.

Setting up the environment

Load the necessary libraries.

library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.0 ──
## ✓ ggplot2 3.3.2     ✓ purrr   0.3.4
## ✓ tibble  3.0.4     ✓ dplyr   1.0.2
## ✓ tidyr   1.1.2     ✓ stringr 1.4.0
## ✓ readr   1.4.0     ✓ forcats 0.5.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(vegan)
## Loading required package: permute
## Loading required package: lattice
## This is vegan 2.5-7
library(ggplot2)
library(cowplot)

Importing data

Import metadata

metadata = read_tsv("year2_sample_mapping_metaphlan.tsv")
## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   .default = col_character(),
##   seq_id_simple = col_double(),
##   fecal_use = col_double(),
##   oral_use = col_double(),
##   mock = col_double(),
##   blank = col_double(),
##   day_number = col_double(),
##   time = col_time(format = ""),
##   extraction_kit_use = col_double(),
##   extraction = col_double(),
##   library_plate = col_double(),
##   id_num = col_double()
## )
## ℹ Use `spec()` for the full column specifications.

Import species-level relative abundance table, transpose, and recalculate relative abundance to sum to 1 instead of 100.

relab = read_tsv("year2_merged_abundance_table_species.txt")
## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   .default = col_double(),
##   body_site = col_character()
## )
## ℹ Use `spec()` for the full column specifications.
relab1 = relab %>% mutate_if(is.numeric, ~./100)
relab1_t = relab1 %>% gather(file, count, 2:164) %>% 
  spread(body_site, count)

Import UniFrac matrices

wuni_fecal = as.dist(
  read.table("year2_merged_abundance_wunifrac_fecal.tsv"))
uuni_fecal = as.dist(
  read.table("year2_merged_abundance_uunifrac_fecal.tsv"))
wuni_oral = as.dist(
  read.table("year2_merged_abundance_wunifrac_oral.tsv"))
uuni_oral = as.dist(
  read.table("year2_merged_abundance_uunifrac_oral.tsv"))

wuni_fecalb = as.dist(
  read.table("year2_merged_abundance_wunifrac_fecalb.tsv"))
uuni_fecalb = as.dist(
  read.table("year2_merged_abundance_uunifrac_fecalb.tsv"))
wuni_oralb = as.dist(
  read.table("year2_merged_abundance_wunifrac_oralb.tsv"))
uuni_oralb = as.dist(
  read.table("year2_merged_abundance_uunifrac_oralb.tsv"))

wuni_fecala = as.dist(
  read.table("year2_merged_abundance_wunifrac_fecala.tsv"))
uuni_fecala = as.dist(
  read.table("year2_merged_abundance_uunifrac_fecala.tsv"))
wuni_orala = as.dist(
  read.table("year2_merged_abundance_wunifrac_orala.tsv"))
uuni_orala = as.dist(
  read.table("year2_merged_abundance_uunifrac_orala.tsv"))

wuni_fecald = as.dist(
  read.table("year2_merged_abundance_wunifrac_fecald.tsv"))
uuni_fecald = as.dist(
  read.table("year2_merged_abundance_uunifrac_fecald.tsv"))

Filter samples

Filter out samples that have no identified species or are controls (blanks and mock communities). Separate oral samples from fecal samples.

relab_fecal = relab1_t %>% left_join(metadata, by = "file") %>% 
  filter(source == "fecal") 
relab_fecal_filtered = relab_fecal %>% 
  mutate(total = rowSums(relab_fecal[,2:496])) %>% 
  filter(total > 0)
relab_oral = relab1_t %>% left_join(metadata, by = "file") %>% 
  filter(source == "oral")

Calculating beta diversity metrics

UniFrac distances are calculated by a MetaPhlAn helper function. Vegan is used to calculate the Bray-Curtis dissimilarity and Jaccard similarity metrics.

brayf = vegdist(relab_fecal_filtered[,2:496], method = "bray")
jaccf = vegdist(relab_fecal_filtered[,2:496], method = "jaccard")
brayo = vegdist(relab_oral[,2:496], method = "bray")
jacco = vegdist(relab_oral[,2:496], method = "jaccard")

PERMANOVAs - All Samples

Differences in the taxonomic structure of the community associated with self-reported ethnicity is assessed using permutational multivariate analysis of variance.

Bray-Curtis

Fecal

adonis2(brayf ~ ethnicity, data = relab_fecal_filtered, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayf ~ ethnicity, data = relab_fecal_filtered, permutations = 4999)
##            Df SumOfSqs      R2      F Pr(>F)    
## ethnicity   1   1.4525 0.06536 7.6225  2e-04 ***
## Residual  109  20.7702 0.93464                  
## Total     110  22.2227 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(brayo ~ ethnicity, data = relab_oral, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayo ~ ethnicity, data = relab_oral, permutations = 4999)
##           Df SumOfSqs     R2      F Pr(>F)   
## ethnicity  1   0.3486 0.0908 3.1957 0.0042 **
## Residual  32   3.4905 0.9092                 
## Total     33   3.8391 1.0000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Bray-Curtis dissimilarity with all samples.

## Run 0 stress 0.1730196 
## Run 1 stress 0.1954393 
## Run 2 stress 0.1890694 
## Run 3 stress 0.1986366 
## Run 4 stress 0.2163515 
## Run 5 stress 0.1730196 
## ... New best solution
## ... Procrustes: rmse 2.113784e-05  max resid 9.989284e-05 
## ... Similar to previous best
## Run 6 stress 0.1857489 
## Run 7 stress 0.1730197 
## ... Procrustes: rmse 9.043795e-05  max resid 0.000386343 
## ... Similar to previous best
## Run 8 stress 0.2305086 
## Run 9 stress 0.2132002 
## Run 10 stress 0.2292344 
## Run 11 stress 0.1878283 
## Run 12 stress 0.1875348 
## Run 13 stress 0.2040644 
## Run 14 stress 0.2066121 
## Run 15 stress 0.1857489 
## Run 16 stress 0.2037089 
## Run 17 stress 0.1931462 
## Run 18 stress 0.17302 
## ... Procrustes: rmse 0.0001923967  max resid 0.0007958677 
## ... Similar to previous best
## Run 19 stress 0.2195464 
## Run 20 stress 0.2084767 
## *** Solution reached

## Run 0 stress 0.1592337 
## Run 1 stress 0.2032347 
## Run 2 stress 0.1592335 
## ... New best solution
## ... Procrustes: rmse 0.000139033  max resid 0.0006605158 
## ... Similar to previous best
## Run 3 stress 0.1754509 
## Run 4 stress 0.1591203 
## ... New best solution
## ... Procrustes: rmse 0.02988077  max resid 0.1460053 
## Run 5 stress 0.1994552 
## Run 6 stress 0.1591203 
## ... Procrustes: rmse 2.196014e-05  max resid 9.723075e-05 
## ... Similar to previous best
## Run 7 stress 0.1665592 
## Run 8 stress 0.1884425 
## Run 9 stress 0.1629265 
## Run 10 stress 0.1660616 
## Run 11 stress 0.1665592 
## Run 12 stress 0.1937613 
## Run 13 stress 0.1632234 
## Run 14 stress 0.1629265 
## Run 15 stress 0.1673417 
## Run 16 stress 0.158507 
## ... New best solution
## ... Procrustes: rmse 0.03767336  max resid 0.1911403 
## Run 17 stress 0.1660616 
## Run 18 stress 0.1766695 
## Run 19 stress 0.158507 
## ... Procrustes: rmse 2.419396e-05  max resid 0.0001101113 
## ... Similar to previous best
## Run 20 stress 0.1591203 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_bray.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(brayf + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(brayf + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(brayo + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Jaccard

Fecal

adonis2(jaccf ~ ethnicity, data = relab_fecal_filtered, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccf ~ ethnicity, data = relab_fecal_filtered, permutations = 4999)
##            Df SumOfSqs      R2      F Pr(>F)    
## ethnicity   1    1.904 0.05969 6.9197  2e-04 ***
## Residual  109   29.992 0.94031                  
## Total     110   31.896 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(jacco ~ ethnicity, data = relab_oral, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jacco ~ ethnicity, data = relab_oral, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)   
## ethnicity  1   0.4911 0.07269 2.5084 0.0042 **
## Residual  32   6.2653 0.92731                 
## Total     33   6.7564 1.00000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Jaccard similarity with all samples.

## Run 0 stress 0.173776 
## Run 1 stress 0.1906485 
## Run 2 stress 0.1919843 
## Run 3 stress 0.201194 
## Run 4 stress 0.1783047 
## Run 5 stress 0.1777442 
## Run 6 stress 0.1957943 
## Run 7 stress 0.2099691 
## Run 8 stress 0.177744 
## Run 9 stress 0.2066552 
## Run 10 stress 0.1730196 
## ... New best solution
## ... Procrustes: rmse 0.03078324  max resid 0.1163285 
## Run 11 stress 0.2185253 
## Run 12 stress 0.1917683 
## Run 13 stress 0.2273878 
## Run 14 stress 0.1880972 
## Run 15 stress 0.2022778 
## Run 16 stress 0.2053309 
## Run 17 stress 0.2215973 
## Run 18 stress 0.2077631 
## Run 19 stress 0.2043171 
## Run 20 stress 0.1928452 
## Run 21 stress 0.2070218 
## Run 22 stress 0.2120901 
## Run 23 stress 0.1838072 
## Run 24 stress 0.1948416 
## Run 25 stress 0.2009606 
## Run 26 stress 0.177744 
## Run 27 stress 0.1765942 
## Run 28 stress 0.2025634 
## Run 29 stress 0.1957946 
## Run 30 stress 0.2033358 
## Run 31 stress 0.1744768 
## Run 32 stress 0.1730196 
## ... New best solution
## ... Procrustes: rmse 1.478532e-05  max resid 7.113282e-05 
## ... Similar to previous best
## *** Solution reached

## Run 0 stress 0.1592336 
## Run 1 stress 0.1820792 
## Run 2 stress 0.1665592 
## Run 3 stress 0.1890765 
## Run 4 stress 0.1593585 
## ... Procrustes: rmse 0.03276102  max resid 0.1485467 
## Run 5 stress 0.1609911 
## Run 6 stress 0.1660763 
## Run 7 stress 0.1861276 
## Run 8 stress 0.1797775 
## Run 9 stress 0.1593585 
## ... Procrustes: rmse 0.03276085  max resid 0.1485464 
## Run 10 stress 0.1609914 
## Run 11 stress 0.1591202 
## ... New best solution
## ... Procrustes: rmse 0.02981504  max resid 0.1456599 
## Run 12 stress 0.1754509 
## Run 13 stress 0.180303 
## Run 14 stress 0.158507 
## ... New best solution
## ... Procrustes: rmse 0.03766872  max resid 0.1911411 
## Run 15 stress 0.1629265 
## Run 16 stress 0.1609932 
## Run 17 stress 0.1665592 
## Run 18 stress 0.1609911 
## Run 19 stress 0.2058922 
## Run 20 stress 0.1819535 
## Run 21 stress 0.1665592 
## Run 22 stress 0.1624179 
## Run 23 stress 0.1927008 
## Run 24 stress 0.1665592 
## Run 25 stress 0.1632234 
## Run 26 stress 0.1665592 
## Run 27 stress 0.1805669 
## Run 28 stress 0.1766695 
## Run 29 stress 0.1593585 
## Run 30 stress 0.162418 
## Run 31 stress 0.1665592 
## Run 32 stress 0.1665592 
## Run 33 stress 0.1806949 
## Run 34 stress 0.1593585 
## Run 35 stress 0.1823253 
## Run 36 stress 0.181956 
## Run 37 stress 0.1629265 
## Run 38 stress 0.1888587 
## Run 39 stress 0.1632234 
## Run 40 stress 0.1819222 
## Run 41 stress 0.1766695 
## Run 42 stress 0.1798959 
## Run 43 stress 0.1660616 
## Run 44 stress 0.201518 
## Run 45 stress 0.1609928 
## Run 46 stress 0.1664842 
## Run 47 stress 0.1629265 
## Run 48 stress 0.1632234 
## Run 49 stress 0.158507 
## ... New best solution
## ... Procrustes: rmse 7.061952e-06  max resid 2.066634e-05 
## ... Similar to previous best
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_jacc.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(jaccf + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(jaccf + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(jacco + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Weighted UniFrac

Fecal

adonis2(wuni_fecal ~ ethnicity, data = relab_fecal_filtered, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_fecal ~ ethnicity, data = relab_fecal_filtered, permutations = 4999)
##            Df SumOfSqs      R2      F Pr(>F)   
## ethnicity   1   0.2429 0.04942 5.6666 0.0012 **
## Residual  109   4.6721 0.95058                 
## Total     110   4.9150 1.00000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(wuni_oral ~ ethnicity, data = relab_oral, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_oral ~ ethnicity, data = relab_oral, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.07415 0.04607 1.5456 0.1736
## Residual  32  1.53533 0.95393              
## Total     33  1.60948 1.00000

NMDS visualization of weighted Unifrac distances with all samples.

## Run 0 stress 0.1128885 
## Run 1 stress 0.1129679 
## ... Procrustes: rmse 0.006469168  max resid 0.0638261 
## Run 2 stress 0.1211153 
## Run 3 stress 0.1129243 
## ... Procrustes: rmse 0.006351517  max resid 0.06336234 
## Run 4 stress 0.1128867 
## ... New best solution
## ... Procrustes: rmse 0.0007263931  max resid 0.005928428 
## ... Similar to previous best
## Run 5 stress 0.1130648 
## ... Procrustes: rmse 0.01034694  max resid 0.08152708 
## Run 6 stress 0.1313283 
## Run 7 stress 0.1378563 
## Run 8 stress 0.1128836 
## ... New best solution
## ... Procrustes: rmse 0.006424424  max resid 0.06385617 
## Run 9 stress 0.1129707 
## ... Procrustes: rmse 0.001786909  max resid 0.01311021 
## Run 10 stress 0.1163974 
## Run 11 stress 0.113063 
## ... Procrustes: rmse 0.007096284  max resid 0.07143052 
## Run 12 stress 0.4130381 
## Run 13 stress 0.1173702 
## Run 14 stress 0.1129006 
## ... Procrustes: rmse 0.006351807  max resid 0.06427751 
## Run 15 stress 0.1421137 
## Run 16 stress 0.1166752 
## Run 17 stress 0.1128898 
## ... Procrustes: rmse 0.0009485802  max resid 0.005774957 
## ... Similar to previous best
## Run 18 stress 0.1128898 
## ... Procrustes: rmse 0.006383938  max resid 0.06402909 
## Run 19 stress 0.1326057 
## Run 20 stress 0.125718 
## *** Solution reached

## Run 0 stress 0.08960119 
## Run 1 stress 0.08968849 
## ... Procrustes: rmse 0.01592553  max resid 0.0774785 
## Run 2 stress 0.1025332 
## Run 3 stress 0.1085128 
## Run 4 stress 0.08968851 
## ... Procrustes: rmse 0.01592196  max resid 0.07746427 
## Run 5 stress 0.103023 
## Run 6 stress 0.1178927 
## Run 7 stress 0.1025332 
## Run 8 stress 0.08968849 
## ... Procrustes: rmse 0.01592304  max resid 0.07746805 
## Run 9 stress 0.0983816 
## Run 10 stress 0.1025332 
## Run 11 stress 0.1157697 
## Run 12 stress 0.1082859 
## Run 13 stress 0.0896012 
## ... Procrustes: rmse 3.168302e-05  max resid 0.0001261021 
## ... Similar to previous best
## Run 14 stress 0.1000692 
## Run 15 stress 0.08968849 
## ... Procrustes: rmse 0.01592722  max resid 0.0774915 
## Run 16 stress 0.1084938 
## Run 17 stress 0.1062022 
## Run 18 stress 0.1172026 
## Run 19 stress 0.1000692 
## Run 20 stress 0.1086048 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_wuni.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(wunif + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(wunif + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(wunio + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Unweighted UniFrac

Fecal

adonis2(uuni_fecal ~ ethnicity, data = relab_fecal_filtered, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_fecal ~ ethnicity, data = relab_fecal_filtered, permutations = 4999)
##            Df SumOfSqs      R2      F Pr(>F)    
## ethnicity   1   0.7732 0.07238 8.5048  2e-04 ***
## Residual  109   9.9100 0.92762                  
## Total     110  10.6832 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(uuni_oral ~ ethnicity, data = relab_oral, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_oral ~ ethnicity, data = relab_oral, permutations = 4999)
##           Df SumOfSqs     R2      F Pr(>F)
## ethnicity  1  0.05327 0.0402 1.3401 0.1916
## Residual  32  1.27193 0.9598              
## Total     33  1.32520 1.0000

NMDS visualization of Bray-Curtis dissimilarity with all samples.

## Run 0 stress 7.26794e-05 
## Run 1 stress 8.89407e-05 
## ... Procrustes: rmse 2.6563e-05  max resid 0.0001140349 
## ... Similar to previous best
## Run 2 stress 9.419245e-05 
## ... Procrustes: rmse 4.891671e-05  max resid 0.0001255239 
## ... Similar to previous best
## Run 3 stress 8.555324e-05 
## ... Procrustes: rmse 2.615283e-05  max resid 0.0001006979 
## ... Similar to previous best
## Run 4 stress 8.908464e-05 
## ... Procrustes: rmse 3.695045e-05  max resid 0.0001250117 
## ... Similar to previous best
## Run 5 stress 8.083382e-05 
## ... Procrustes: rmse 3.973176e-05  max resid 0.0001347788 
## ... Similar to previous best
## Run 6 stress 9.079105e-05 
## ... Procrustes: rmse 3.511466e-05  max resid 0.0001306709 
## ... Similar to previous best
## Run 7 stress 9.147368e-05 
## ... Procrustes: rmse 2.977313e-05  max resid 0.0001073563 
## ... Similar to previous best
## Run 8 stress 7.159914e-05 
## ... New best solution
## ... Procrustes: rmse 2.972201e-05  max resid 0.0001033728 
## ... Similar to previous best
## Run 9 stress 8.509091e-05 
## ... Procrustes: rmse 3.107911e-05  max resid 9.833484e-05 
## ... Similar to previous best
## Run 10 stress 7.854869e-05 
## ... Procrustes: rmse 2.14154e-05  max resid 7.716503e-05 
## ... Similar to previous best
## Run 11 stress 8.615579e-05 
## ... Procrustes: rmse 2.458636e-05  max resid 8.14055e-05 
## ... Similar to previous best
## Run 12 stress 7.366728e-05 
## ... Procrustes: rmse 1.979365e-05  max resid 9.986138e-05 
## ... Similar to previous best
## Run 13 stress 9.662353e-05 
## ... Procrustes: rmse 2.424657e-05  max resid 6.917625e-05 
## ... Similar to previous best
## Run 14 stress 9.077723e-05 
## ... Procrustes: rmse 4.125896e-05  max resid 0.000110526 
## ... Similar to previous best
## Run 15 stress 6.001949e-05 
## ... New best solution
## ... Procrustes: rmse 2.02629e-05  max resid 6.462539e-05 
## ... Similar to previous best
## Run 16 stress 7.398229e-05 
## ... Procrustes: rmse 1.990127e-05  max resid 6.056812e-05 
## ... Similar to previous best
## Run 17 stress 7.391904e-05 
## ... Procrustes: rmse 2.237521e-05  max resid 6.279285e-05 
## ... Similar to previous best
## Run 18 stress 9.026151e-05 
## ... Procrustes: rmse 3.250961e-05  max resid 7.56077e-05 
## ... Similar to previous best
## Run 19 stress 9.061626e-05 
## ... Procrustes: rmse 2.611029e-05  max resid 0.0001082571 
## ... Similar to previous best
## Run 20 stress 8.268287e-05 
## ... Procrustes: rmse 2.122476e-05  max resid 6.927945e-05 
## ... Similar to previous best
## *** Solution reached
## Warning in metaMDS(uuni_fecal, k = 2, trymax = 500): stress is (nearly) zero:
## you may have insufficient data
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Run 0 stress 0.1791554 
## Run 1 stress 0.1819031 
## Run 2 stress 0.1879203 
## Run 3 stress 0.1834742 
## Run 4 stress 0.1833687 
## Run 5 stress 0.1830513 
## Run 6 stress 0.1825586 
## Run 7 stress 0.2191732 
## Run 8 stress 0.1933048 
## Run 9 stress 0.1806554 
## Run 10 stress 0.193551 
## Run 11 stress 0.1815046 
## Run 12 stress 0.1791029 
## ... New best solution
## ... Procrustes: rmse 0.007986632  max resid 0.03645454 
## Run 13 stress 0.1791029 
## ... New best solution
## ... Procrustes: rmse 3.331413e-05  max resid 9.393246e-05 
## ... Similar to previous best
## Run 14 stress 0.1795051 
## ... Procrustes: rmse 0.01244262  max resid 0.05107072 
## Run 15 stress 0.1791029 
## ... Procrustes: rmse 5.581876e-05  max resid 0.0001815511 
## ... Similar to previous best
## Run 16 stress 0.1909943 
## Run 17 stress 0.1873035 
## Run 18 stress 0.1792178 
## ... Procrustes: rmse 0.01941215  max resid 0.053807 
## Run 19 stress 0.1814713 
## Run 20 stress 0.1791554 
## ... Procrustes: rmse 0.007993025  max resid 0.03651193 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_uuni.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(uunif + theme(legend.box.margin = margin(0, 0, 0, 12)))
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
col1 = plot_grid(uunif + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
col2 = plot_grid(uunio + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

PERMANOVAs - Before samples

Filtering samples

Filter samples from before the diet intervention and recalculate Bray-Curtis and Jaccard metrics

relab_fecal_filtered_before = relab_fecal_filtered %>% 
  filter(Study_period == "Before")
relab_oral_before = relab_oral %>% 
  filter(Study_period == "Before")

brayfb = vegdist(relab_fecal_filtered_before[,2:496], method = "bray")
jaccfb = vegdist(relab_fecal_filtered_before[,2:496], method = "jaccard")
brayob = vegdist(relab_oral_before[,2:496], method = "bray")
jaccob = vegdist(relab_oral_before[,2:496], method = "jaccard")

Bray-Curtis

Fecal

adonis2(brayfb ~ ethnicity, data = relab_fecal_filtered_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayfb ~ ethnicity, data = relab_fecal_filtered_before, permutations = 4999)
##           Df SumOfSqs      R2     F Pr(>F)
## ethnicity  1   0.2744 0.06875 1.255 0.2314
## Residual  17   3.7171 0.93125             
## Total     18   3.9915 1.00000

Oral

adonis2(brayob ~ ethnicity, data = relab_oral_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayob ~ ethnicity, data = relab_oral_before, permutations = 4999)
##           Df SumOfSqs      R2     F Pr(>F)
## ethnicity  1  0.18873 0.09118 1.505 0.1208
## Residual  15  1.88112 0.90882             
## Total     16  2.06985 1.00000

NMDS visualization of Bray-Curtis dissimilarity with before samples.

## Run 0 stress 0.1284112 
## Run 1 stress 0.1712267 
## Run 2 stress 0.1444672 
## Run 3 stress 0.1284112 
## ... Procrustes: rmse 5.107342e-06  max resid 1.072505e-05 
## ... Similar to previous best
## Run 4 stress 0.1612495 
## Run 5 stress 0.1284112 
## ... Procrustes: rmse 5.270567e-06  max resid 1.099011e-05 
## ... Similar to previous best
## Run 6 stress 0.1284112 
## ... New best solution
## ... Procrustes: rmse 1.868258e-06  max resid 3.636777e-06 
## ... Similar to previous best
## Run 7 stress 0.1617329 
## Run 8 stress 0.1842236 
## Run 9 stress 0.2014128 
## Run 10 stress 0.1444672 
## Run 11 stress 0.1506507 
## Run 12 stress 0.1712268 
## Run 13 stress 0.1712267 
## Run 14 stress 0.1444672 
## Run 15 stress 0.1792601 
## Run 16 stress 0.158196 
## Run 17 stress 0.1618869 
## Run 18 stress 0.1842236 
## Run 19 stress 0.1525329 
## Run 20 stress 0.158196 
## *** Solution reached

## Run 0 stress 0.1253357 
## Run 1 stress 0.1481009 
## Run 2 stress 0.1453195 
## Run 3 stress 0.1364623 
## Run 4 stress 0.1585964 
## Run 5 stress 0.1591175 
## Run 6 stress 0.1453195 
## Run 7 stress 0.1585963 
## Run 8 stress 0.1498923 
## Run 9 stress 0.1453195 
## Run 10 stress 0.1591175 
## Run 11 stress 0.1481009 
## Run 12 stress 0.1403222 
## Run 13 stress 0.1300515 
## Run 14 stress 0.1405455 
## Run 15 stress 0.1364625 
## Run 16 stress 0.3689216 
## Run 17 stress 0.1253357 
## ... Procrustes: rmse 0.0001814997  max resid 0.0005916947 
## ... Similar to previous best
## Run 18 stress 0.1253358 
## ... Procrustes: rmse 0.0002713796  max resid 0.0008824368 
## ... Similar to previous best
## Run 19 stress 0.1686967 
## Run 20 stress 0.1467916 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_bray_before.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(brayfb + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(brayfb + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(brayob + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Jaccard

Fecal

adonis2(jaccfb ~ ethnicity, data = relab_fecal_filtered_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccfb ~ ethnicity, data = relab_fecal_filtered_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1   0.4122 0.07421 1.3628 0.1328
## Residual  17   5.1416 0.92579              
## Total     18   5.5537 1.00000

Oral

adonis2(jaccob ~ ethnicity, data = relab_oral_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccob ~ ethnicity, data = relab_oral_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1   0.2815 0.07968 1.2988 0.1338
## Residual  15   3.2508 0.92032              
## Total     16   3.5323 1.00000

NMDS visualization of Jaccard similarity with before samples.

## Run 0 stress 0.1284112 
## Run 1 stress 0.1444672 
## Run 2 stress 0.161887 
## Run 3 stress 0.1518819 
## Run 4 stress 0.1284112 
## ... Procrustes: rmse 4.827122e-06  max resid 9.395989e-06 
## ... Similar to previous best
## Run 5 stress 0.1518819 
## Run 6 stress 0.1284112 
## ... Procrustes: rmse 2.297949e-06  max resid 4.477037e-06 
## ... Similar to previous best
## Run 7 stress 0.162817 
## Run 8 stress 0.1284112 
## ... Procrustes: rmse 3.074062e-06  max resid 6.583554e-06 
## ... Similar to previous best
## Run 9 stress 0.1581961 
## Run 10 stress 0.1842236 
## Run 11 stress 0.1617329 
## Run 12 stress 0.2626449 
## Run 13 stress 0.1284112 
## ... Procrustes: rmse 3.799604e-06  max resid 7.560866e-06 
## ... Similar to previous best
## Run 14 stress 0.1612495 
## Run 15 stress 0.1284112 
## ... Procrustes: rmse 2.600561e-06  max resid 5.256494e-06 
## ... Similar to previous best
## Run 16 stress 0.1791334 
## Run 17 stress 0.1792602 
## Run 18 stress 0.1712269 
## Run 19 stress 0.1525329 
## Run 20 stress 0.1284112 
## ... Procrustes: rmse 4.904981e-06  max resid 1.056895e-05 
## ... Similar to previous best
## *** Solution reached

## Run 0 stress 0.1253358 
## Run 1 stress 0.1591172 
## Run 2 stress 0.1364623 
## Run 3 stress 0.1455176 
## Run 4 stress 0.1295583 
## Run 5 stress 0.1689227 
## Run 6 stress 0.1253357 
## ... New best solution
## ... Procrustes: rmse 0.0001331197  max resid 0.0004425188 
## ... Similar to previous best
## Run 7 stress 0.1256166 
## ... Procrustes: rmse 0.02144201  max resid 0.0619013 
## Run 8 stress 0.1461135 
## Run 9 stress 0.1295583 
## Run 10 stress 0.1405455 
## Run 11 stress 0.1384923 
## Run 12 stress 0.1453195 
## Run 13 stress 0.3633597 
## Run 14 stress 0.1403222 
## Run 15 stress 0.1364623 
## Run 16 stress 0.1300761 
## Run 17 stress 0.1591176 
## Run 18 stress 0.1484361 
## Run 19 stress 0.1364624 
## Run 20 stress 0.1467916 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_jacc_before.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(jaccfb + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(jaccfb + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(jaccob + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Weighted UniFrac

Fecal

adonis2(wuni_fecalb ~ ethnicity, data = relab_fecal_filtered_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_fecalb ~ ethnicity, data = relab_fecal_filtered_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.02942 0.03414 0.6009 0.6444
## Residual  17  0.83220 0.96586              
## Total     18  0.86162 1.00000

Oral

adonis2(wuni_oralb ~ ethnicity, data = relab_oral_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_oralb ~ ethnicity, data = relab_oral_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.04241 0.04763 0.7501 0.4686
## Residual  15  0.84808 0.95237              
## Total     16  0.89049 1.00000

NMDS visualization of weighted Unifrac distances with before samples.

## Run 0 stress 0.07991057 
## Run 1 stress 0.2321205 
## Run 2 stress 0.08273949 
## Run 3 stress 0.07999755 
## ... Procrustes: rmse 0.009288069  max resid 0.02907652 
## Run 4 stress 0.08403423 
## Run 5 stress 0.08148164 
## Run 6 stress 0.1424107 
## Run 7 stress 0.08148173 
## Run 8 stress 0.1475728 
## Run 9 stress 0.08273945 
## Run 10 stress 0.07991036 
## ... New best solution
## ... Procrustes: rmse 0.0004527546  max resid 0.001336798 
## ... Similar to previous best
## Run 11 stress 0.08148166 
## Run 12 stress 0.08273955 
## Run 13 stress 0.08273944 
## Run 14 stress 0.07999764 
## ... Procrustes: rmse 0.009378223  max resid 0.02906354 
## Run 15 stress 0.08403422 
## Run 16 stress 0.07991039 
## ... Procrustes: rmse 0.0003635338  max resid 0.00100614 
## ... Similar to previous best
## Run 17 stress 0.08273945 
## Run 18 stress 0.08273949 
## Run 19 stress 0.1453388 
## Run 20 stress 0.08403423 
## *** Solution reached

## Run 0 stress 0.06894997 
## Run 1 stress 0.06894997 
## ... New best solution
## ... Procrustes: rmse 8.094821e-06  max resid 2.596814e-05 
## ... Similar to previous best
## Run 2 stress 0.07415714 
## Run 3 stress 0.06966296 
## Run 4 stress 0.07415714 
## Run 5 stress 0.08144912 
## Run 6 stress 0.06966296 
## Run 7 stress 0.06966296 
## Run 8 stress 0.0761536 
## Run 9 stress 0.08571785 
## Run 10 stress 0.07047078 
## Run 11 stress 0.06894997 
## ... New best solution
## ... Procrustes: rmse 5.478821e-06  max resid 1.561733e-05 
## ... Similar to previous best
## Run 12 stress 0.0704708 
## Run 13 stress 0.07643534 
## Run 14 stress 0.07415715 
## Run 15 stress 0.07415714 
## Run 16 stress 0.07643538 
## Run 17 stress 0.06894998 
## ... Procrustes: rmse 3.032653e-05  max resid 9.395695e-05 
## ... Similar to previous best
## Run 18 stress 0.285 
## Run 19 stress 0.06894998 
## ... Procrustes: rmse 1.473e-05  max resid 3.59209e-05 
## ... Similar to previous best
## Run 20 stress 0.06894997 
## ... Procrustes: rmse 6.664391e-06  max resid 1.804576e-05 
## ... Similar to previous best
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_wuni_before.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(wunifb + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(wunifb + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(wuniob + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Unweighted UniFrac

Fecal

adonis2(uuni_fecalb ~ ethnicity, data = relab_fecal_filtered_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_fecalb ~ ethnicity, data = relab_fecal_filtered_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.18751 0.08793 1.6389 0.0182 *
## Residual  17  1.94502 0.91207                
## Total     18  2.13253 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(uuni_oralb ~ ethnicity, data = relab_oral_before, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_oralb ~ ethnicity, data = relab_oral_before, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.03204 0.04827 0.7608 0.7076
## Residual  15  0.63162 0.95173              
## Total     16  0.66365 1.00000

NMDS visualization of unweighted UniFrac distances with before samples.

## Run 0 stress 9.874988e-05 
## Run 1 stress 9.961836e-05 
## ... Procrustes: rmse 4.474348e-05  max resid 8.783092e-05 
## ... Similar to previous best
## Run 2 stress 7.274105e-05 
## ... New best solution
## ... Procrustes: rmse 0.0001072195  max resid 0.00024281 
## ... Similar to previous best
## Run 3 stress 8.935065e-05 
## ... Procrustes: rmse 9.023071e-05  max resid 0.000223894 
## ... Similar to previous best
## Run 4 stress 9.359378e-05 
## ... Procrustes: rmse 0.0001046243  max resid 0.0002039388 
## ... Similar to previous best
## Run 5 stress 9.251685e-05 
## ... Procrustes: rmse 0.0001186332  max resid 0.00025548 
## ... Similar to previous best
## Run 6 stress 9.371519e-05 
## ... Procrustes: rmse 0.0001033606  max resid 0.0002122299 
## ... Similar to previous best
## Run 7 stress 9.45746e-05 
## ... Procrustes: rmse 0.0001015814  max resid 0.0002177958 
## ... Similar to previous best
## Run 8 stress 8.808076e-05 
## ... Procrustes: rmse 8.321261e-05  max resid 0.000185678 
## ... Similar to previous best
## Run 9 stress 8.642016e-05 
## ... Procrustes: rmse 5.583815e-05  max resid 0.0001016477 
## ... Similar to previous best
## Run 10 stress 9.755127e-05 
## ... Procrustes: rmse 9.972745e-05  max resid 0.0002585689 
## ... Similar to previous best
## Run 11 stress 9.511959e-05 
## ... Procrustes: rmse 9.932362e-05  max resid 0.0002001329 
## ... Similar to previous best
## Run 12 stress 8.969811e-05 
## ... Procrustes: rmse 0.0001205089  max resid 0.0002574019 
## ... Similar to previous best
## Run 13 stress 8.796747e-05 
## ... Procrustes: rmse 8.56072e-05  max resid 0.0002052298 
## ... Similar to previous best
## Run 14 stress 9.693815e-05 
## ... Procrustes: rmse 0.0001012259  max resid 0.0002085296 
## ... Similar to previous best
## Run 15 stress 8.442425e-05 
## ... Procrustes: rmse 0.0001122393  max resid 0.0002247436 
## ... Similar to previous best
## Run 16 stress 9.699694e-05 
## ... Procrustes: rmse 0.0001184309  max resid 0.0002528589 
## ... Similar to previous best
## Run 17 stress 9.119991e-05 
## ... Procrustes: rmse 9.086308e-05  max resid 0.0002544837 
## ... Similar to previous best
## Run 18 stress 9.839301e-05 
## ... Procrustes: rmse 0.0001046333  max resid 0.0002329734 
## ... Similar to previous best
## Run 19 stress 9.280435e-05 
## ... Procrustes: rmse 0.0001143868  max resid 0.0002475763 
## ... Similar to previous best
## Run 20 stress 8.419121e-05 
## ... Procrustes: rmse 6.22314e-05  max resid 0.0001601429 
## ... Similar to previous best
## *** Solution reached
## Warning in metaMDS(uuni_fecalb, k = 2, trymax = 500): stress is (nearly) zero:
## you may have insufficient data
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure

## Run 0 stress 0.1694422 
## Run 1 stress 0.1703134 
## Run 2 stress 0.1767621 
## Run 3 stress 0.1710156 
## Run 4 stress 0.1694423 
## ... Procrustes: rmse 6.34574e-05  max resid 0.0001699156 
## ... Similar to previous best
## Run 5 stress 0.1697768 
## ... Procrustes: rmse 0.03161416  max resid 0.08143387 
## Run 6 stress 0.1768018 
## Run 7 stress 0.1715766 
## Run 8 stress 0.1710154 
## Run 9 stress 0.203838 
## Run 10 stress 0.1767621 
## Run 11 stress 0.2174096 
## Run 12 stress 0.2133062 
## Run 13 stress 0.1767621 
## Run 14 stress 0.1694423 
## ... Procrustes: rmse 8.271933e-05  max resid 0.0002249821 
## ... Similar to previous best
## Run 15 stress 0.1703132 
## Run 16 stress 0.1694423 
## ... Procrustes: rmse 3.127789e-05  max resid 8.070865e-05 
## ... Similar to previous best
## Run 17 stress 0.1767623 
## Run 18 stress 0.1710162 
## Run 19 stress 0.2038383 
## Run 20 stress 0.1710157 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_uuni_before.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(uunifb + theme(legend.box.margin = margin(0, 0, 0, 12)))
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
col1 = plot_grid(uunifb + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
## Warning in MASS::cov.trob(data[, vars]): Probable convergence failure
col2 = plot_grid(uuniob + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

PERMANOVAs - During samples

Filter samples

Filter samples from during the diet intervention and recalculate Bray-Curtis and Jaccard metrics

relab_fecal_filtered_during = relab_fecal_filtered %>% 
  filter(Study_period == "During")

brayfd = vegdist(relab_fecal_filtered_during[,2:496], method = "bray")
jaccfd = vegdist(relab_fecal_filtered_during[,2:496], method = "jaccard")

Bray-Curtis

Fecal

adonis2(brayfd ~ ethnicity, data = relab_fecal_filtered_during, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayfd ~ ethnicity, data = relab_fecal_filtered_during, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)    
## ethnicity  1   0.7899 0.06703 4.2387  2e-04 ***
## Residual  59  10.9957 0.93297                  
## Total     60  11.7856 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Bray-Curtis dissimilarity with during samples.

## Run 0 stress 0.1664788 
## Run 1 stress 0.1900921 
## Run 2 stress 0.2004255 
## Run 3 stress 0.1900921 
## Run 4 stress 0.2004255 
## Run 5 stress 0.2002677 
## Run 6 stress 0.2007708 
## Run 7 stress 0.1900406 
## Run 8 stress 0.1752755 
## Run 9 stress 0.1664617 
## ... New best solution
## ... Procrustes: rmse 0.001066307  max resid 0.006181484 
## ... Similar to previous best
## Run 10 stress 0.1664617 
## ... Procrustes: rmse 3.476369e-05  max resid 0.000119953 
## ... Similar to previous best
## Run 11 stress 0.2248416 
## Run 12 stress 0.1664617 
## ... Procrustes: rmse 8.390135e-05  max resid 0.0002723608 
## ... Similar to previous best
## Run 13 stress 0.2187408 
## Run 14 stress 0.2195501 
## Run 15 stress 0.1664617 
## ... Procrustes: rmse 8.367231e-06  max resid 2.661781e-05 
## ... Similar to previous best
## Run 16 stress 0.208741 
## Run 17 stress 0.1664788 
## ... Procrustes: rmse 0.001068306  max resid 0.006189377 
## ... Similar to previous best
## Run 18 stress 0.2003321 
## Run 19 stress 0.1900406 
## Run 20 stress 0.2266991 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_bray_during.tif", res = 150, width = 11, 
     height = 7, units = "in")
brayfd
dev.off()
## quartz_off_screen 
##                 2

Jaccard

Fecal

adonis2(jaccfd ~ ethnicity, data = relab_fecal_filtered_during, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccfd ~ ethnicity, data = relab_fecal_filtered_during, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)    
## ethnicity  1   1.0741 0.06259 3.9392  4e-04 ***
## Residual  59  16.0879 0.93741                  
## Total     60  17.1620 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Jaccard similarity with during samples.

## Run 0 stress 0.1664788 
## Run 1 stress 0.172702 
## Run 2 stress 0.223669 
## Run 3 stress 0.1664788 
## ... New best solution
## ... Procrustes: rmse 4.915646e-05  max resid 0.0001542918 
## ... Similar to previous best
## Run 4 stress 0.1965621 
## Run 5 stress 0.19664 
## Run 6 stress 0.1664617 
## ... New best solution
## ... Procrustes: rmse 0.001068985  max resid 0.006190718 
## ... Similar to previous best
## Run 7 stress 0.2008191 
## Run 8 stress 0.1900394 
## Run 9 stress 0.1664617 
## ... Procrustes: rmse 6.804985e-05  max resid 0.0002170154 
## ... Similar to previous best
## Run 10 stress 0.2007708 
## Run 11 stress 0.2004655 
## Run 12 stress 0.1784168 
## Run 13 stress 0.2039778 
## Run 14 stress 0.1724188 
## Run 15 stress 0.1664617 
## ... Procrustes: rmse 7.791776e-06  max resid 2.582068e-05 
## ... Similar to previous best
## Run 16 stress 0.1727032 
## Run 17 stress 0.1900973 
## Run 18 stress 0.1664617 
## ... Procrustes: rmse 1.324225e-05  max resid 4.279012e-05 
## ... Similar to previous best
## Run 19 stress 0.2039778 
## Run 20 stress 0.1998369 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_jacc_during.tif", res = 150, width = 11, 
     height = 7, units = "in")
jaccfd
dev.off()
## quartz_off_screen 
##                 2

Weighted UniFrac

Fecal

adonis2(wuni_fecald ~ ethnicity, data = relab_fecal_filtered_during, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_fecald ~ ethnicity, data = relab_fecal_filtered_during, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.15963 0.06033 3.7881 0.0134 *
## Residual  59  2.48636 0.93967                
## Total     60  2.64599 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of weighted Unifrac distances with during samples.

## Run 0 stress 0.1065766 
## Run 1 stress 0.1066173 
## ... Procrustes: rmse 0.004441149  max resid 0.02908068 
## Run 2 stress 0.121613 
## Run 3 stress 0.1170402 
## Run 4 stress 0.1347412 
## Run 5 stress 0.1066765 
## ... Procrustes: rmse 0.0117531  max resid 0.08209334 
## Run 6 stress 0.1328284 
## Run 7 stress 0.110359 
## Run 8 stress 0.1157495 
## Run 9 stress 0.110029 
## Run 10 stress 0.1421085 
## Run 11 stress 0.1171166 
## Run 12 stress 0.1065766 
## ... New best solution
## ... Procrustes: rmse 5.599846e-05  max resid 0.0002443808 
## ... Similar to previous best
## Run 13 stress 0.1159533 
## Run 14 stress 0.1066172 
## ... Procrustes: rmse 0.004443165  max resid 0.02906563 
## Run 15 stress 0.1381071 
## Run 16 stress 0.1378626 
## Run 17 stress 0.1157515 
## Run 18 stress 0.138338 
## Run 19 stress 0.1182267 
## Run 20 stress 0.1125892 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_wuni_during.tif", res = 150, width = 11, 
     height = 7, units = "in")
wunif
dev.off()
## quartz_off_screen 
##                 2

Unweighted UniFrac

Fecal

adonis2(uuni_fecald ~ ethnicity, data = relab_fecal_filtered_during, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_fecald ~ ethnicity, data = relab_fecal_filtered_during, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)    
## ethnicity  1   0.4436 0.07779 4.9764  2e-04 ***
## Residual  59   5.2595 0.92221                  
## Total     60   5.7031 1.00000                  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of unweighted Unifrac distance with during samples.

## Run 0 stress 0.1793941 
## Run 1 stress 0.1815779 
## Run 2 stress 0.1793928 
## ... New best solution
## ... Procrustes: rmse 0.001474469  max resid 0.008378042 
## ... Similar to previous best
## Run 3 stress 0.1799881 
## Run 4 stress 0.1800008 
## Run 5 stress 0.1794325 
## ... Procrustes: rmse 0.005388824  max resid 0.02821562 
## Run 6 stress 0.1794325 
## ... Procrustes: rmse 0.00536919  max resid 0.02815032 
## Run 7 stress 0.2121887 
## Run 8 stress 0.1866832 
## Run 9 stress 0.1799881 
## Run 10 stress 0.179654 
## ... Procrustes: rmse 0.01012475  max resid 0.03493015 
## Run 11 stress 0.1847463 
## Run 12 stress 0.1813386 
## Run 13 stress 0.1793933 
## ... Procrustes: rmse 0.001043957  max resid 0.005966139 
## ... Similar to previous best
## Run 14 stress 0.1796595 
## ... Procrustes: rmse 0.01082075  max resid 0.03796503 
## Run 15 stress 0.2119015 
## Run 16 stress 0.2120047 
## Run 17 stress 0.2184909 
## Run 18 stress 0.1800008 
## Run 19 stress 0.2171491 
## Run 20 stress 0.1794952 
## ... Procrustes: rmse 0.007252083  max resid 0.03211737 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_uuni_during.tif", res = 150, width = 18, 
     height = 7, units = "in")
uunifd
dev.off()
## quartz_off_screen 
##                 2

PERMANOVAs - After samples

Filter samples

Filter after samples and recalculate Bray-Curtis and Jaccard metrics

relab_fecal_filtered_after = relab_fecal_filtered %>% 
  filter(Study_period == "After")
relab_oral_after = relab_oral %>% 
  filter(Study_period == "After")

brayfa = vegdist(relab_fecal_filtered_after[,2:496], method = "bray")
jaccfa = vegdist(relab_fecal_filtered_after[,2:496], method = "jaccard")
brayoa = vegdist(relab_oral_after[,2:496], method = "bray")
jaccoa = vegdist(relab_oral_after[,2:496], method = "jaccard")

Bray-Curtis

Fecal

adonis2(brayfa ~ ethnicity, data = relab_fecal_filtered_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayfa ~ ethnicity, data = relab_fecal_filtered_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1   0.4762 0.07579 2.3782 0.0208 *
## Residual  29   5.8065 0.92421                
## Total     30   6.2827 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(brayoa ~ ethnicity, data = relab_oral_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = brayoa ~ ethnicity, data = relab_oral_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.24789 0.14292 2.5014  0.018 *
## Residual  15  1.48655 0.85708                
## Total     16  1.73444 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Bray-Curtis dissimilarity with after samples.

## Run 0 stress 0.1590347 
## Run 1 stress 0.1632431 
## Run 2 stress 0.1711856 
## Run 3 stress 0.1912912 
## Run 4 stress 0.1590551 
## ... Procrustes: rmse 0.004813471  max resid 0.01859784 
## Run 5 stress 0.1909684 
## Run 6 stress 0.1632431 
## Run 7 stress 0.1590347 
## ... New best solution
## ... Procrustes: rmse 2.462918e-05  max resid 7.695649e-05 
## ... Similar to previous best
## Run 8 stress 0.2009304 
## Run 9 stress 0.1916395 
## Run 10 stress 0.1913422 
## Run 11 stress 0.1590347 
## ... New best solution
## ... Procrustes: rmse 1.255126e-05  max resid 3.775422e-05 
## ... Similar to previous best
## Run 12 stress 0.1590347 
## ... Procrustes: rmse 1.23226e-05  max resid 5.312786e-05 
## ... Similar to previous best
## Run 13 stress 0.2053616 
## Run 14 stress 0.1885132 
## Run 15 stress 0.1904057 
## Run 16 stress 0.1802135 
## Run 17 stress 0.1656964 
## Run 18 stress 0.2359138 
## Run 19 stress 0.1590347 
## ... Procrustes: rmse 1.633959e-05  max resid 5.295247e-05 
## ... Similar to previous best
## Run 20 stress 0.1632431 
## *** Solution reached

## Run 0 stress 0.1367903 
## Run 1 stress 0.1367906 
## ... Procrustes: rmse 0.0001540849  max resid 0.0003625897 
## ... Similar to previous best
## Run 2 stress 0.1345672 
## ... New best solution
## ... Procrustes: rmse 0.03762177  max resid 0.1126176 
## Run 3 stress 0.1345672 
## ... New best solution
## ... Procrustes: rmse 7.113183e-06  max resid 1.54181e-05 
## ... Similar to previous best
## Run 4 stress 0.1367906 
## Run 5 stress 0.1573084 
## Run 6 stress 0.1367904 
## Run 7 stress 0.1367904 
## Run 8 stress 0.1520581 
## Run 9 stress 0.1367902 
## Run 10 stress 0.1345672 
## ... New best solution
## ... Procrustes: rmse 2.539426e-06  max resid 4.973324e-06 
## ... Similar to previous best
## Run 11 stress 0.1563232 
## Run 12 stress 0.2845416 
## Run 13 stress 0.1563232 
## Run 14 stress 0.1504065 
## Run 15 stress 0.1367904 
## Run 16 stress 0.1345672 
## ... Procrustes: rmse 1.34276e-05  max resid 2.875464e-05 
## ... Similar to previous best
## Run 17 stress 0.1504065 
## Run 18 stress 0.1367902 
## Run 19 stress 0.1769399 
## Run 20 stress 0.1577923 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_bray_after.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(brayfa + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(brayfa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(brayoa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Jaccard

Fecal

adonis2(jaccfa ~ ethnicity, data = relab_fecal_filtered_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccfa ~ ethnicity, data = relab_fecal_filtered_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1   0.6328 0.07096 2.2149 0.0144 *
## Residual  29   8.2848 0.92904                
## Total     30   8.9175 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(jaccoa ~ ethnicity, data = relab_oral_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = jaccoa ~ ethnicity, data = relab_oral_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.36814 0.11723 1.9919 0.0168 *
## Residual  15  2.77222 0.88277                
## Total     16  3.14036 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

NMDS visualization of Jaccard similarity with after samples.

## Run 0 stress 0.1590347 
## Run 1 stress 0.1776623 
## Run 2 stress 0.3892089 
## Run 3 stress 0.1966504 
## Run 4 stress 0.2027067 
## Run 5 stress 0.159055 
## ... Procrustes: rmse 0.004810536  max resid 0.01859397 
## Run 6 stress 0.1590554 
## ... Procrustes: rmse 0.0048899  max resid 0.0188668 
## Run 7 stress 0.1590547 
## ... Procrustes: rmse 0.004675963  max resid 0.01804715 
## Run 8 stress 0.1773823 
## Run 9 stress 0.2218271 
## Run 10 stress 0.2014925 
## Run 11 stress 0.1909684 
## Run 12 stress 0.1590347 
## ... Procrustes: rmse 2.903439e-05  max resid 8.579623e-05 
## ... Similar to previous best
## Run 13 stress 0.1590347 
## ... Procrustes: rmse 2.174493e-05  max resid 8.522461e-05 
## ... Similar to previous best
## Run 14 stress 0.1590348 
## ... Procrustes: rmse 0.0001101951  max resid 0.0002671833 
## ... Similar to previous best
## Run 15 stress 0.1866588 
## Run 16 stress 0.1878211 
## Run 17 stress 0.2009304 
## Run 18 stress 0.1978081 
## Run 19 stress 0.1590347 
## ... Procrustes: rmse 1.189067e-05  max resid 2.539754e-05 
## ... Similar to previous best
## Run 20 stress 0.1632431 
## *** Solution reached

## Run 0 stress 0.1367902 
## Run 1 stress 0.1345672 
## ... New best solution
## ... Procrustes: rmse 0.03748147  max resid 0.1125011 
## Run 2 stress 0.1520582 
## Run 3 stress 0.1520583 
## Run 4 stress 0.1367902 
## Run 5 stress 0.1345672 
## ... Procrustes: rmse 2.018597e-05  max resid 4.552282e-05 
## ... Similar to previous best
## Run 6 stress 0.1769413 
## Run 7 stress 0.1504065 
## Run 8 stress 0.1367902 
## Run 9 stress 0.1504065 
## Run 10 stress 0.1563232 
## Run 11 stress 0.1345672 
## ... Procrustes: rmse 1.494667e-05  max resid 3.577065e-05 
## ... Similar to previous best
## Run 12 stress 0.1345672 
## ... Procrustes: rmse 2.801263e-05  max resid 6.683795e-05 
## ... Similar to previous best
## Run 13 stress 0.1504066 
## Run 14 stress 0.1367903 
## Run 15 stress 0.1345672 
## ... Procrustes: rmse 3.525083e-05  max resid 9.642475e-05 
## ... Similar to previous best
## Run 16 stress 0.1345672 
## ... Procrustes: rmse 2.232049e-05  max resid 4.431062e-05 
## ... Similar to previous best
## Run 17 stress 0.1563232 
## Run 18 stress 0.1367905 
## Run 19 stress 0.1637332 
## Run 20 stress 0.1367905 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_jacc_after.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(jaccfa + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(jaccfa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(jaccoa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Weighted UniFrac

Fecal

adonis2(wuni_fecala ~ ethnicity, data = relab_fecal_filtered_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_fecala ~ ethnicity, data = relab_fecal_filtered_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)  
## ethnicity  1  0.09394 0.07012 2.1868 0.0766 .
## Residual  29  1.24577 0.92988                
## Total     30  1.33971 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(wuni_orala ~ ethnicity, data = relab_oral_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = wuni_orala ~ ethnicity, data = relab_oral_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.07966 0.11412 1.9323 0.1274
## Residual  15  0.61840 0.88588              
## Total     16  0.69806 1.00000

NMDS visualization of weighted Unifrac distances with after samples.

## Run 0 stress 0.103874 
## Run 1 stress 0.1259496 
## Run 2 stress 0.1030343 
## ... New best solution
## ... Procrustes: rmse 0.03047841  max resid 0.1470687 
## Run 3 stress 0.1079067 
## Run 4 stress 0.1030343 
## ... Procrustes: rmse 7.058277e-05  max resid 0.0002370727 
## ... Similar to previous best
## Run 5 stress 0.1030343 
## ... New best solution
## ... Procrustes: rmse 1.996087e-05  max resid 6.624552e-05 
## ... Similar to previous best
## Run 6 stress 0.1038765 
## Run 7 stress 0.1147664 
## Run 8 stress 0.1079067 
## Run 9 stress 0.1079067 
## Run 10 stress 0.112468 
## Run 11 stress 0.1038764 
## Run 12 stress 0.1079067 
## Run 13 stress 0.3906131 
## Run 14 stress 0.1038693 
## Run 15 stress 0.1030343 
## ... Procrustes: rmse 1.588753e-05  max resid 5.537614e-05 
## ... Similar to previous best
## Run 16 stress 0.1283419 
## Run 17 stress 0.1030343 
## ... Procrustes: rmse 6.011413e-05  max resid 0.0002865742 
## ... Similar to previous best
## Run 18 stress 0.1030343 
## ... New best solution
## ... Procrustes: rmse 4.031004e-06  max resid 1.078038e-05 
## ... Similar to previous best
## Run 19 stress 0.1030343 
## ... New best solution
## ... Procrustes: rmse 2.101056e-05  max resid 9.684372e-05 
## ... Similar to previous best
## Run 20 stress 0.1038692 
## *** Solution reached

## Run 0 stress 0.05926746 
## Run 1 stress 0.09412422 
## Run 2 stress 0.07447768 
## Run 3 stress 0.08351463 
## Run 4 stress 0.0866882 
## Run 5 stress 0.09687134 
## Run 6 stress 0.07191275 
## Run 7 stress 0.05926746 
## ... Procrustes: rmse 1.264535e-06  max resid 3.067386e-06 
## ... Similar to previous best
## Run 8 stress 0.07401827 
## Run 9 stress 0.0866882 
## Run 10 stress 0.08521648 
## Run 11 stress 0.07447768 
## Run 12 stress 0.08726034 
## Run 13 stress 0.08605518 
## Run 14 stress 0.07191278 
## Run 15 stress 0.07447768 
## Run 16 stress 0.05926746 
## ... Procrustes: rmse 2.169967e-05  max resid 6.844801e-05 
## ... Similar to previous best
## Run 17 stress 0.05926746 
## ... New best solution
## ... Procrustes: rmse 6.442793e-06  max resid 1.983388e-05 
## ... Similar to previous best
## Run 18 stress 0.05926746 
## ... Procrustes: rmse 5.890497e-06  max resid 1.692085e-05 
## ... Similar to previous best
## Run 19 stress 0.09687134 
## Run 20 stress 0.07191276 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_wuni_after.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(wunifa + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(wunifa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(wunioa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2

Unweighted UniFrac

Fecal

adonis2(uuni_fecala ~ ethnicity, data = relab_fecal_filtered_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_fecala ~ ethnicity, data = relab_fecal_filtered_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)   
## ethnicity  1  0.20279 0.07361 2.3043 0.0052 **
## Residual  29  2.55215 0.92639                 
## Total     30  2.75494 1.00000                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Oral

adonis2(uuni_orala ~ ethnicity, data = relab_oral_after, 
        permutations = 4999)
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 4999
## 
## adonis2(formula = uuni_orala ~ ethnicity, data = relab_oral_after, permutations = 4999)
##           Df SumOfSqs      R2      F Pr(>F)
## ethnicity  1  0.04386 0.06782 1.0914 0.3476
## Residual  15  0.60285 0.93218              
## Total     16  0.64671 1.00000

NMDS visualization of unweighted UniFrac distances with after samples.

## Run 0 stress 0.1657305 
## Run 1 stress 0.1658233 
## ... Procrustes: rmse 0.00537419  max resid 0.02533288 
## Run 2 stress 0.1833738 
## Run 3 stress 0.2137625 
## Run 4 stress 0.1789747 
## Run 5 stress 0.1781267 
## Run 6 stress 0.1686699 
## Run 7 stress 0.1678211 
## Run 8 stress 0.1681404 
## Run 9 stress 0.1824032 
## Run 10 stress 0.1657192 
## ... New best solution
## ... Procrustes: rmse 0.001665992  max resid 0.005852336 
## ... Similar to previous best
## Run 11 stress 0.1681404 
## Run 12 stress 0.1675407 
## Run 13 stress 0.1657194 
## ... Procrustes: rmse 0.0001596051  max resid 0.0007250493 
## ... Similar to previous best
## Run 14 stress 0.1686699 
## Run 15 stress 0.1835957 
## Run 16 stress 0.1658235 
## ... Procrustes: rmse 0.00589096  max resid 0.02558163 
## Run 17 stress 0.1686699 
## Run 18 stress 0.1789746 
## Run 19 stress 0.1675406 
## Run 20 stress 0.1657189 
## ... New best solution
## ... Procrustes: rmse 0.0003874098  max resid 0.001768522 
## ... Similar to previous best
## *** Solution reached

## Run 0 stress 0.1364074 
## Run 1 stress 0.1694069 
## Run 2 stress 0.1442196 
## Run 3 stress 0.1442203 
## Run 4 stress 0.1304387 
## ... New best solution
## ... Procrustes: rmse 0.04532598  max resid 0.1276811 
## Run 5 stress 0.1380133 
## Run 6 stress 0.1304387 
## ... New best solution
## ... Procrustes: rmse 3.720389e-05  max resid 0.0001027526 
## ... Similar to previous best
## Run 7 stress 0.1543537 
## Run 8 stress 0.136732 
## Run 9 stress 0.1304387 
## ... Procrustes: rmse 1.198968e-05  max resid 3.044183e-05 
## ... Similar to previous best
## Run 10 stress 0.1670458 
## Run 11 stress 0.1439974 
## Run 12 stress 0.1442198 
## Run 13 stress 0.1367313 
## Run 14 stress 0.1384723 
## Run 15 stress 0.1304386 
## ... New best solution
## ... Procrustes: rmse 8.04969e-05  max resid 0.0002287135 
## ... Similar to previous best
## Run 16 stress 0.1442208 
## Run 17 stress 0.1442203 
## Run 18 stress 0.1679828 
## Run 19 stress 0.1445522 
## Run 20 stress 0.1571375 
## *** Solution reached

tiff(file = "figures/nmds_plot_taxa_uuni_after.tif", res = 150, width = 18, 
     height = 7, units = "in")
legend1 = get_legend(uunifa + theme(legend.box.margin = margin(0, 0, 0, 12)))
col1 = plot_grid(uunifa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('A'), 
                 label_size = 20)
col2 = plot_grid(uunioa + theme(legend.position = "none"),
                 nrow = 1, ncol = 1, labels = c('B'), 
                 label_size = 20) 
plot = plot_grid(col1, col2, legend1, nrow = 1, ncol = 3, 
          rel_widths = c(1.5, 1.5, 0.75), align = "h", axis = "t")
plot
dev.off()
## quartz_off_screen 
##                 2